<<<<<<< HEAD
"""
黄河小浪底调水调沙问题
"""
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
from scipy.interpolate import CubicSpline
from scipy.optimize import curve_fit


print(__doc__)
plt.style.use(['grid', 'muted'])

# 数据处理
data = pd.read_csv('data/XiaoLangDi.csv', sep=',')
data['time'] = [3600*(12*i-4) for i in range(1, len(data)+1, 1)]
data['sand discharge'] = data['water-flow rate'] * data['sediment']

t0 = data['time'].values
y0 = data['sand discharge'].values

# 利用插值估计任意时刻的排沙量
cs = CubicSpline(t0, y0, bc_type='natural')
t1 = np.linspace(t0.min(), t0.max(), 50)
y1 = cs(t1)

# 利用拟合研究排沙量和水流量的关系
func = lambda x, a, b, c: a*np.power(x, 2) + b*x + c
s0 = data['sediment'].values
w0 = data['water-flow rate'].values

popt, _ = curve_fit(func, w0, s0)
a, b, c = popt
s1 = func(w0, a, b, c)

# 可视化
plt.figure(figsize=(8, 6))

plt.subplot(1, 2, 1)
plt.scatter(t0, y0, c='b', label='origin')
plt.plot(t1, y1, 'r--', label='interpolation')
plt.title('XiaoLangDi Sand Discharge')
plt.xlabel('t')
plt.ylabel('y')
plt.legend(loc="upper right")

plt.subplot(1, 2, 2)
plt.plot(w0, s0, '*', label='original values')
plt.plot(w0, s1, 'r--', label='curve_fit values')
plt.title('S-W Curve')
plt.xlabel('w')
plt.ylabel('s')
plt.legend(loc="upper right")

plt.show()
=======
"""
黄河小浪底调水调沙问题
"""
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
from scipy.interpolate import CubicSpline
from scipy.optimize import curve_fit


print(__doc__)
plt.style.use(['grid', 'muted'])

# 数据处理
data = pd.read_csv('data/XiaoLangDi.csv', sep=',')
data['time'] = [3600*(12*i-4) for i in range(1, len(data)+1, 1)]
data['sand discharge'] = data['water-flow rate'] * data['sediment']

t0 = data['time'].values
y0 = data['sand discharge'].values

# 利用插值估计任意时刻的排沙量
cs = CubicSpline(t0, y0, bc_type='natural')
t1 = np.linspace(t0.min(), t0.max(), 50)
y1 = cs(t1)

# 利用拟合研究排沙量和水流量的关系
func = lambda x, a, b, c: a*np.power(x, 2) + b*x + c
s0 = data['sediment'].values
w0 = data['water-flow rate'].values

popt, _ = curve_fit(func, w0, s0)
a, b, c = popt
s1 = func(w0, a, b, c)

# 可视化
plt.figure(figsize=(8, 6))

plt.subplot(1, 2, 1)
plt.scatter(t0, y0, c='b', label='origin')
plt.plot(t1, y1, 'r--', label='interpolation')
plt.title('XiaoLangDi Sand Discharge')
plt.xlabel('t')
plt.ylabel('y')
plt.legend(loc="upper right")

plt.subplot(1, 2, 2)
plt.plot(w0, s0, '*', label='original values')
plt.plot(w0, s1, 'r--', label='curve_fit values')
plt.title('S-W Curve')
plt.xlabel('w')
plt.ylabel('s')
plt.legend(loc="upper right")

plt.show()
>>>>>>> a66c8eec2c3bbe955d7da215f43ffffda9c7b6b5
